COVID-19 Progression

library(tidyverse)
library(magrittr)
library(lubridate)
library(scales)
library(matrixStats)
library(ggrepel)
library(broom)
library(glue)
library(jsonlite)
library(rvest)
library(RCurl)
library(pander)
library(DT)
library(plotly)
library(cowplot)
library(QuantTools)
panderOptions("big.mark", ",")
panderOptions("table.split.table", Inf)
panderOptions("table.style", "rmarkdown")
panderOptions("missing", "")
theme_set(theme_bw())

# Handle updates between 12am & 9am
dt <- Sys.Date()
if (as.numeric(format(Sys.time(), "%H")) <= 9){
  dt <- Sys.Date() - 1
}

auStates <- c(
  ACT = "Australian Capital Territory",
  QLD = "Queensland",
  NSW = "New South Wales",
  VIC = "Victoria",
  SA = "South Australia",
  WA = "Western Australia",
  NT = "Northern Territory",
  TAS = "Tasmania"
)

Disclaimer: This very simple report was prepared by a bioinformatician with no experience in epidemiology or virology, and as such should be treated simply as an alternate viewpoint on the data, which I was simply unable to find elsewhere. Many other people exist with much greater expertise on this subject. However, I do hope this provides a useful perspective which is able to add constructively to the wider discussion. In addition, it should be noted that this is very much focussed on Australian data.

Data Sources

confirmed <- url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_covid19_confirmed_global.csv") %>%
    read_csv() %>%
    pivot_longer(
        cols = ends_with("20"),
        names_to = "date",
        values_to = "confirmed"
    )  %>%
    mutate(
        date = str_replace_all(
            date, "(.+)/(.+)/(.+)", "20\\3-\\1-\\2"
        ) %>%
            ymd()
    ) %>%
    dplyr::rename(
        Country = `Country/Region`
    ) %>%
    dplyr::mutate(
        Country = case_when(
            `Province/State` == "Hubei" ~ "China (Hubei)",
            `Province/State` == "Hong Kong" ~ "Hong Kong",
            grepl("China", Country) & !`Province/State` %in% c("Hubei", "Hong Kong") ~ "China (Other)",
            Country == "Korea, South" ~ "South Korea",
            Country == "Congo (Kinshasa)" ~ "DR Congo",
            Country == "Congo (Brazzaville)" ~ "Congo",
            !grepl("China", Country) ~ Country
        )
    ) %>%
    dplyr::filter(
        !is.na(confirmed)
    ) %>%
  dplyr::select(-Lat, -Long)

guardData <- fromJSON("https://interactive.guim.co.uk/docsdata/1q5gdePANXci8enuiS4oHUJxcxC13d6bjMRSicakychE.json")
altAU <- guardData$sheets$updates %>% 
  as_tibble() %>% 
  mutate(
    `Province/State` = auStates[State], 
    Country = "Australia", 
    Date = parse_date_time(Date, orders = "%d/%m/%Y") %>% ymd()
  ) %>% 
  dplyr::select(`Province/State`, Country, date = Date, confirmed = `Cumulative case count`) %>% 
  mutate(confirmed = as.numeric(confirmed)) %>%
  arrange(`Province/State`, date) %>%
  tidyr::complete(`Province/State`, Country, date) %>%
  group_by(`Province/State`) %>%
  fill(confirmed)  %>%
  dplyr::filter(!is.na(confirmed))  %>%
  ungroup()
# Add the data, with the Guardian API being preferentialy selected
confirmed %<>%
  mutate(source = "JHU") %>%
  bind_rows(
    mutate(altAU, source = "API")
  ) %>%
  arrange(
    date, `Province/State`, source
  ) %>%
  distinct(`Province/State`, Country, date, .keep_all = TRUE) %>%
  dplyr::select(-source)

latestAU <- list()

nswRecov <- nswTests <- NA_real_
nswUrl <- "https://www.health.nsw.gov.au/_layouts/feed.aspx?xsl=1&web=/news&page=4ac47e14-04a9-4016-b501-65a23280e841&wp=baabf81e-a904-44f1-8d59-5f6d56519965&pageurl=/news/Pages/rss-nsw-health.aspx" %>%
  read_xml() %>% 
  xml_find_all("//item") %>% 
  .[[1]] %>%
  as.character() %>% 
  str_split("\n") %>% 
  .[[1]] %>% 
  str_subset("https") %>% 
  str_extract("https[^<]+")
nswTable <- nswUrl %>% 
  read_html() %>% 
  html_nodes("body") %>% 
  xml_find_all("//table") %>% 
  .[[1]] %>% 
  html_table() %>%
  set_colnames(
     str_replace_all(colnames(.), "\\u200b", "")
  ) %>%
  mutate_all(str_replace_all, pattern = "\\u200b", replacement = "") %>%
  mutate(
    Count = str_remove_all(Count, "[\\*,]") %>% as.numeric,
    Cases = str_replace_all(Cases, " +", " ")
  )
nswTests <- nswTable %>%
  dplyr::filter(
    str_detect(Cases, "Total tests carried out")
  ) %>%
  .[["Count"]]
nswRecov <- nswTable %>%
  dplyr::filter(str_detect(Cases, "recovered")) %>%
  .[["Count"]] 
latestAU$NSW <- tibble(
  State = "New South Wales",
  date = dt,
  confirmed = nswTable %>% 
    dplyr::filter(str_detect(Cases, "Confirmed cases")) %>%
    .[["Count"]],
  deaths = nswTable %>% 
    dplyr::filter(str_detect(Cases, "Deaths")) %>%
    .[["Count"]],
  recovered = nswRecov,
  tested = nswTests
)

qldRecov <- qldTests <- NA_real_
qldUrl <- "https://www.qld.gov.au/health/conditions/health-alerts/coronavirus-covid-19/current-status/statistics"
qldTests <- qldUrl %>%
  read_html() %>%
  html_nodes("body") %>%
  xml_find_all("//div[contains(@class, 'tested')]") %>%
  html_text() %>%
  str_extract("[0-9,]+") %>%
  str_remove_all(",") %>%
  as.numeric()
qldTable <- qldUrl %>%
  read_html() %>%
  html_nodes("body") %>%
  xml_find_all("//table[@id = 'QLD_Cases']") %>%
  html_table() %>%
  .[[1]] 
qldRecov <- qldTable %>%
  dplyr::filter(Cases == "Recovered") %>%
  .[["Total"]] %>%
  str_remove_all(",") %>%
  as.numeric()
latestAU$QLD <- tibble(
  State = "Queensland",
  date = dt,
  confirmed = qldTable %>% 
    dplyr::filter(str_detect(Cases, "Number of (cases|confirmed)")) %>%
    .[["Total"]] %>% 
    str_remove_all(",") %>%
    str_extract("[0-9]+") %>% 
    as.numeric(),
  deaths = qldTable %>% 
    dplyr::filter(str_detect(Cases, "Deaths")) %>%
    .[["Total"]] %>% 
    str_remove_all(",") %>%
    str_extract("[0-9]+") %>% 
    as.numeric(),
  recovered = qldRecov,
  tested = qldTests
)

vicRecov <- vicTests <- NA_real_
vicUrl <- c(
  "https://www.dhhs.vic.gov.au/coronavirus-update-victoria-{format(dt, '%d')}-{format.Date(dt, '%B')}-{year(dt)}",
  "https://www.dhhs.vic.gov.au/coronavirus-update-victoria-{format(dt - 1, '%d')}-{format.Date(dt - 1, '%B')}-{year(dt-1)}"
  ) %>% 
  vapply(glue, character(1)) %>%
  str_to_lower() %>%
  vapply(url.exists, logical(1)) %>%
  which() %>%
  names() %>%
  .[[1]] 
vicTxt <- vicUrl %>%
  read_html() %>%
  html_nodes("body") %>% 
  html_text() %>% 
  str_split("\n") %>% 
  .[[1]] %>% 
  str_trim() %>% 
  .[. != ""] %>%
  str_subset("(test|recovered|total|died)") %>% 
  stringi::stri_trans_general("latin-ascii")
vicTests <- vicTxt %>%
  str_subset(".*[0-9]+.*tests") %>%
  .[[1]] %>%
  str_replace_all(
    ".* ([0-9,]+)[^0-9]test(ed|s|ing).+",
    "\\1"
  ) %>%
  str_remove_all(",") %>%
  as.numeric()
vicRecov <- vicTxt %>%
  str_subset("[0-9]+.+recovered") %>%
  str_trim() %>%
  str_replace_all(".*[^0-9,]([0-9,]+)[^0-9]*people have recovered.*", "\\1") %>%
  str_remove(",") %>%
  as.numeric() 
latestAU$VIC <- tibble(
  State = "Victoria",
  date = dt,
  confirmed = vicTxt %>%
    str_subset("(coronavirus|total).+cases") %>% 
    str_subset("Victoria") %>%
    str_subset("regional", negate = TRUE) %>%
    str_replace_all(".+[^0-9,]+([0-9,]+)[^0-9,].+", "\\1") %>% 
    str_remove_all(",") %>%
    as.numeric() %>%
    unique(),
  deaths = vicTxt %>% 
    str_subset("died") %>% 
    str_replace_all(".+[^0-9]([0-9,]+)[^0-9]*people have died.+", "\\1") %>%
    str_remove_all(",") %>% 
    as.numeric() %>%
    unique(),
  recovered = vicRecov,
  tested = vicTests
)

waTests <- waRecov <- NA_real_
waMedia <- "https://ww2.health.wa.gov.au/News/Media-releases-listing-page" %>%
  read_html() %>% 
  html_nodes("body") %>% 
  xml_find_all("//a[contains(@title, '19 update')]") %>%
  .[[1]] %>%
  html_attr("href")
waUrl <- glue("https://ww2.health.wa.gov.au{waMedia}")
waTxt <- waUrl %>%
  read_html() %>%
  html_nodes("body") %>%
  html_text() %>%
  str_split("\n") %>%
  .[[1]] %>%
  str_trim() %>%
  .[. != ""]
waTests <- waTxt %>%
  str_subset("test(s|ed).*performed") %>%
  str_replace_all(".+ ([0-9,]+) COVID.+", "\\1") %>%
  str_remove_all(",") %>%
  as.numeric()
waRecov <- waTxt %>%
  str_subset("recovered") %>%
  str_subset("[0-9]..") %>%
  str_remove_all("COVID.*19") %>% 
  str_replace_all("[^0-9]*([0-9]*)[^0-9]*", "\\1") %>%
  as.numeric()
latestAU$WA <- tibble(
  State = "Western Australia",
  date = dt,
  confirmed =waTxt %>% 
    stringi::stri_trans_general("latin-ascii") %>% 
    str_subset("State's") %>% 
    str_remove_all("COVID.*19") %>% 
    str_replace_all("[^0-9]*([0-9]+)[^0-9]*", "\\1") %>% 
    as.numeric(),
  deaths = 9, # Needs to be done manually again
  recovered = waRecov,
  tested = waTests
)

saRecov <- saTests <- NA_real_
saUrl <- "https://www.covid-19.sa.gov.au/home/dashboard"
saTxt <- saUrl  %>%
  read_html() %>%
  xml_find_all("//div[contains(@class, 'accordion__content')]") %>%
  .[[1]] %>%
  html_text() %>%
  str_split("\r\n") %>%
  .[[1]] %>%
  str_trim() %>%
  .[. != ""]
## The dashboard is unreliable so use SA Health press releases as an additional source
## These don't appear until late afternoon, so this is a pain in the arse quite frankly
altSaUrl <- c(
  "https://www.sahealth.sa.gov.au/wps/wcm/connect/public+content/sa+health+internet/about+us/news+and+media/all+media+releases/covid-19+update+{day(dt)}+{format(dt, '%B')}+{year(dt)}",
  "https://www.sahealth.sa.gov.au/wps/wcm/connect/public+content/sa+health+internet/about+us/news+and+media/all+media+releases/covid-19+update+{day(dt - 1)}+{format(dt - 1, '%B')}+{year(dt)}"
  ) %>% 
  vapply(glue, character(1)) %>% 
  str_to_lower()
validSa <- altSaUrl %>%
  lapply(read_html) %>%
  lapply(html_node, css = "body") %>%
  vapply(function(x){!is.na(x)}, logical(1)) %>%
  which() %>%
  .[[1]]
altSaTxt <- altSaUrl[validSa] %>%
  read_html() %>%
  html_node("body") %>%
  xml_find_all("//section[@class = 'main-content']") %>%
  html_text() %>%
  str_split("\n") %>%
  .[[1]] %>%
  str_trim() %>%
  .[. != ""]
saTests <- c(
  saTxt %>%
    str_replace_all(".+(conducted|undertaken).*[^0-9,]([0-9, ]+)[^0-9,]+tests.+", "\\2") %>%
    str_remove_all("[, ]") %>%
    as.numeric(),
  altSaTxt %>% 
    str_subset("tests") %>% 
    str_replace_all(".+ ([0-9,]+) test.+", "\\1") %>% 
    str_remove_all(",") %>% 
    as.numeric()
) %>%
  max(na.rm = TRUE)
saTable <- "https://www.covid-19.sa.gov.au/" %>% 
  read_html() %>% 
  html_node("body") %>% 
  xml_find_all("//div[@class = 'focus-boxes']") %>%
  html_text() %>% 
  str_split("\r\n") %>% 
  .[[1]] %>% 
  str_trim() %>% 
  .[.!=""] %>% 
  .[seq(length(.), 1, by = -1)] %>% 
  .[!str_detect(., "(Data current|View the|See more)")] %>%
  matrix(ncol = 2, byrow = TRUE) %>% 
  set_colnames(c("Cases", "Total")) %>% 
  as.data.frame(stringsAsFactors = FALSE) %>% 
  mutate(
    Total = as.numeric(Total),
    Cases = str_remove_all(Cases, "in South Australia")
  )
saRecov <- c(
  saTable %>%
    dplyr::filter(str_detect(Cases, "recovered")) %>%
    .[["Total"]],
  altSaTxt %>% 
    str_subset("(cleared|recovered)") %>% 
    str_replace_all(".+ ([0-9]+) people have.+", "\\1") %>% 
    as.numeric()
) %>%
  max(na.rm = TRUE)
latestAU$SA <- tibble(
  State = "South Australia",
  date = dt,
  confirmed = saTable %>% 
    dplyr::filter(str_detect(Cases, "Total cases")) %>% 
    .[["Total"]],
  deaths = 4, # Only able to be added manually. Thanks SA Health
  recovered = saRecov,
  tested = saTests
)

actRecov <- actTests <- NA_real_
actUrl <- "https://www.covid19.act.gov.au/home"
actTable <- actUrl %>%
  read_html() %>%
  html_nodes("body") %>%
  xml_find_all(
    "//div[contains(@class, 'spf-article-card--tabular')]"
  ) %>% 
  html_text() %>% 
  str_split("\r\n") %>% 
  .[[4]] %>% 
  str_trim() %>% 
  setdiff(y = "") %>% 
  str_remove_all("[,*]") %>%
  matrix(ncol = 2, byrow = TRUE) %>%
  set_colnames(c("Cases", "Total")) %>%
  as.data.frame(stringsAsFactors = FALSE) %>%
  mutate(Total = as.numeric(Total))
actTests <- actTable %>% 
  dplyr::filter(!str_detect(Cases, "recovered")) %>%
  .[["Total"]] %>%
  sum()
actRecov <- actTable %>%
  dplyr::filter(str_detect(Cases, "recov")) %>%
  .[["Total"]]
latestAU$ACT <- tibble(
  State = "Australian Capital Territory",
  date = dt,
  confirmed = actTable %>% 
    dplyr::filter(str_detect(Cases, "Confirmed cases")) %>% 
    .[["Total"]],
  deaths = 3, # Only able to be added manually. Thanks SA Health
  recovered = actRecov,
  tested = actTests
)

tasUrl <- "https://www.coronavirus.tas.gov.au/facts/cases-and-testing-updates"
tasTables <- tasUrl %>%
  read_html() %>%
  html_nodes("body") %>%
  xml_find_all("//table") %>%
  html_table() %>%
  .[1:2] %>%
  lapply(mutate, Number = str_remove_all(Number, "[^0-9]")) %>%
  lapply(mutate, Number  = as.numeric(Number))
latestAU$TAS <- tibble(
  State = "Tasmania",
  date = dt,
  confirmed = dplyr::filter(
     tasTables[[1]], 
     str_detect(`Cases in Tasmania`, "Total cases")
     )$Number,
  deaths = dplyr::filter(
     tasTables[[1]], 
     str_detect(`Cases in Tasmania`, "Deaths")
     )$Number,
  recovered = dplyr::filter(
    tasTables[[1]],
    `Cases in Tasmania` == "Recovered"
    )$Number,
  tested = dplyr::filter(
    tasTables[[2]],
    `Laboratory tests` == "Total laboratory tests"
    )$Number
)

ntRecov <- ntTests <- NA_real_
ntUrl <- "https://coronavirus.nt.gov.au/"
ntTxt <- ntUrl %>%
  read_html() %>%
  html_nodes("body") %>%
  xml_find_all(
    "//div[contains(@class, 'header-widget')]"
  ) %>%
  html_text() %>%
    str_split("\r\n") %>%
    .[[1]] %>%
    str_trim() %>%
    .[. != ""]
ntTests <- ntTxt %>%
  str_subset("test") %>%
  str_extract("([0-9,]+)") %>%
  str_remove_all(",") %>%
  as.numeric()
ntRecov <- ntTxt %>%
  str_subset("recovered") %>%
  str_extract("[0-9,]+") %>%
  str_remove_all(",") %>%
  as.numeric()
latestAU$NT <- tibble(
  State = "Northern Territory",
  date = dt,
  confirmed = 30,# ntTxt %>%
    # str_subset("confirmed") %>% 
    # str_extract("[0-9]+") %>% 
    # as.numeric(),
  deaths = 0,
  recovered = ntRecov,
  tested = ntTests
)

latestAU %<>% 
  bind_rows() %>%
  mutate(Country = "Australia")

latestAU %>%
  dplyr::select(
    `Province/State` = State, Country, date, recovered
  ) %>%
  write_tsv(
    here::here(glue("recovered/recovered_{dt}.tsv"))
  )

confirmed %<>%
  bind_rows(
    dplyr::select(latestAU, any_of(colnames(.)), `Province/State` = State)
  ) %>%
  arrange(date) %>%
  mutate(f = paste(`Province/State`, Country, sep = "_")) %>%
  split(f = .$f) %>%
  lapply(mutate, confirmed = cummax(confirmed)) %>%
  bind_rows() %>%
  dplyr::select(-f) %>%
  arrange(Country, `Province/State`, date, desc(confirmed)) %>%
  dplyr::distinct(`Province/State`, Country, date, .keep_all = TRUE)

auRecTS <- list.files(here::here("recovered"), pattern = "rec.*tsv",full.names = TRUE) %>%
  lapply(read_tsv) %>%
  bind_rows()

recovered <- url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_covid19_recovered_global.csv") %>%
  read_csv() %>%
  pivot_longer(
    cols = ends_with("20"),
    names_to = "date",
    values_to = "recovered"
  )  %>%
  mutate(
    date = str_replace_all(
      date, "(.+)/(.+)/(.+)", "20\\3-\\1-\\2"
    ) %>%
      ymd()
  ) %>%
  dplyr::rename(
    Country = `Country/Region`
  ) %>%
  dplyr::mutate(
    Country = case_when(
      `Province/State` == "Hubei" ~ "China (Hubei)",
      `Province/State` == "Hong Kong" ~ "Hong Kong",
      grepl("China", Country) & !`Province/State` %in% c("Hubei", "Hong Kong") ~ "China (Other)",
      grepl("Korea, South", Country) ~ "South Korea",
      Country == "Congo (Kinshasa)" ~ "DR Congo",
      Country == "Congo (Brazzaville)" ~ "Congo",
      !grepl("China", Country) ~ Country
    )
  ) %>%
  dplyr::select(-Lat, -Long) %>%
  bind_rows(auRecTS) %>%
  group_by(
    `Province/State`, Country, date
  ) %>%
  summarise_at(vars(recovered), max)

deaths <- url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_covid19_deaths_global.csv") %>%
  read_csv() %>%
  pivot_longer(
    cols = ends_with("20"),
    names_to = "date",
    values_to = "deaths"
  )  %>%
  mutate(
    date = str_replace_all(
      date, "(.+)/(.+)/(.+)", "20\\3-\\1-\\2"
    ) %>%
      ymd()
  ) %>%
  dplyr::rename(
    Country = `Country/Region`
  ) %>%
  dplyr::mutate(
    Country = case_when(
      `Province/State` == "Hubei" ~ "China (Hubei)",
      `Province/State` == "Hong Kong" ~ "Hong Kong",
      grepl("China", Country) & !`Province/State` %in% c("Hubei", "Hong Kong") ~ "China (Other)",
      Country == "Korea, South" ~ "South Korea",
      Country == "Congo (Kinshasa)" ~ "DR Congo",
      Country == "Congo (Brazzaville)" ~ "Congo",      
      !grepl("China", Country) ~ Country
    )
  ) %>%
  dplyr::select(-Lat, -Long) %>%
  bind_rows(
    dplyr::select(latestAU, any_of(colnames(.)), `Province/State` = State)
  ) %>%
  arrange(date) %>%
  mutate(f = paste(`Province/State`, Country, sep = "_")) %>%
  split(f = .$f) %>%
  lapply(mutate, deaths = cummax(deaths)) %>%
  bind_rows() %>%
  dplyr::select(-f) %>%
  dplyr::distinct(`Province/State`, Country, date, .keep_all = TRUE)

wikiPops <- read_tsv("wikiPops.tsv") %>%
  mutate_at(c("Region", "Continent"), as.factor)

alwaysShow <- c("Australia", "New Zealand", "US", "United Kingdom", "Taiwan*", "Singapore", "South Korea", "China (Other)", "China (Hubei)", "Hong Kong", "Japan", "Ireland", "Russia", "Brazil")

Data for confirmed cases, recoveries and fatalities was primarily sourced from Johns Hopkins University, using the datasets provided at https://github.com/CSSEGISandData/COVID-19. JHU data is now updated every 24 hours at approximately 3:30(UTC), which is about 1:00PM in Adelaide. As such no accurate, daily updates for international data can be produced until after that time. Importantly, dates associated with confirmed cases from this data source, may differ from dates associated with confirmed cases from Australian sources. For example, cases reported in the morning in Australia may be assigned to the previous day in US data sources.

Live hourly updates for Australia are available from https://covid-19-au.github.io/ for those who would like an up to the minute breakdown of confirmed cases. Numbers used for generation of this page are updated periodically throughout the day using values provided by individual states at NSW, QLD, VIC, WA, SA, TAS, ACT and the NT. Additional Australian data was obtained from The Guardian. Whilst dates from the Guardian are likely to be reported using the UK as a reference, this is still closer to Australian time zones than USA-based data.

Population sizes were obtained from 2019 UN Estimates. Given the disparity of infection within China, China was broken into Hubei Province and the rest of China, with Hong Kong and Taiwan already being considered as separate countries in all datasets. Population estimates for Hubei Province were taken from the 2018 estimates given by Statista.com and this is likely to be a very slight underestimate.

Confirmed cases of COVID-19 as provided by the Chinese Government have been discussed elsewhere as unusual, and data appears potentially unreliable. In this analysis, discussions regarding accurate Chinese reporting are not considered further and data is simply presented as supplied by JHU.

However, all countries are likely to contain many unreported cases given the incomplete testing regimes in place for most countries. Similarly, reporting in many countries may have features that cause concerns regarding data integrity and this makes comparison across countries difficult. Information on recovered cases has been difficult to accurately obtain due inconsistent methods for considering a case as recovered, and lack of reporting for these cases in many jurisdictions. In Australia, some states (e.g. NSW & QLD) have only begun releasing these figures in mid-April, whilst other states such as Victoria were releasing these numbers from mid-March.

International Data

startingPoint <- 4
minPop <- 2e6

For this section, most data is presented relative to population size. Growth in infection rates is only shown after the point at which the cumulative confirmed infection rate breached 4 confirmed cases / million. This equates to about 101 confirmed cases within Australia, and is broadly comparable to the “Days since passing 100 confirmed cases” commonly shown elsewhere.

Regional Cumulative Infection Rates

Before exploring individual countries, a regional perspective may be helpful. Importantly, only Oceania, Eastern Asia, Western Europe and Southern Europe have successfully controlled the infection rate. Most other regions are still experiencing exponential growth. The time period for the most recent doubling of global cases was 35 days. At the time of report preparation, the total number of confirmed cases, including all countries for which data has been made available, is 6,736,474.

ggplotly(
  confirmed %>% 
    inner_join(wikiPops) %>%
    group_by(Continent, Region, Country, Population, date) %>% 
    summarise(confirmed = sum(confirmed)) %>% 
    ungroup() %>%
    group_by(Continent, Region, date) %>% 
    summarise(rate = round(1e6*sum(confirmed)/sum(Population), 2)) %>% 
    dplyr::filter(date < dt, date > "2020-02-18", rate > 1) %>%
    ggplot(aes(date, rate, colour = Region)) +  
    geom_line() +
    scale_y_log10(label = comma_format(1)) + 
    facet_wrap(~Continent, ncol = 2) +
    labs(x = "Date", y = "Cumulative Infection Rate (Confirmed Cases / Million)")
)

Current regional cumulative infection rates. Lines are shown from the point that local cases exceed 1 confirmed case / million

Cumulative Infection Rates

Table of Most Impacted Countries

Confirmed cases in this table are effectively the cumulative, confirmed incidence rate. Recovered patients and those who have passed away are still included in these numbers.

confirmed %>%
  group_by(Country, date) %>%
  summarise(confirmed = sum(confirmed)) %>%
  ungroup() %>%
  group_by(Country) %>%
  dplyr::filter(
    date == max(date),
  ) %>%
  ungroup() %>%
  inner_join(wikiPops) %>%
  mutate(
    rate = 1e6*confirmed / Population,
    occurrence = Population / confirmed,
    Population = round(Population, -3) / 1e6
  ) %>%
  dplyr::filter(rate >= 1) %>%
  arrange(desc(rate)) %>%
  rename_at(vars(everything()), str_to_title) %>%
  dplyr::select(
    Continent, Region, Country,
    Date, Confirmed, 
    `Population (millions)` = Population, 
    `Rate (Cases per Million)` = Rate,
    Occurrence
  ) %>%
  datatable(
    options = list(
      pageLength = 25, 
      autoWidth = TRUE,
      searchCols = list(
        NULL, NULL, NULL, NULL, NULL,
        list(
          search = glue(
            '{minPop/1e6 + 0.001} ... {max(wikiPops$Population/1e6)}'
          )
        ),
        NULL
      )
    ),
    filter = 'top',
    class = "stripe",
    rownames = FALSE,
    caption = paste(
      "The most impacted countries studied here and shown as a proportion of total population.",
      "The initial filter is set so that only countries with a population greater than", comma(minPop), "are shown.",
      "All fields are searchable and sortable.",
      "To filter numeric columns, either use the slider or enter the values in the form 'min ... max'.",
      "To filter text columns, partial matching used in a case-insensitive manner.",
      "Populations have been rounded to the nearest thousand to make reading values easier.",
      "'Rate' represents the latest confirmed infection rate as cumulative cases per million people, whilst 'Occurrence' represents the number of people expected before one case is found, assuming an even distribution amongst the population.",
      "In other words, one in every 'Occurrence' people within the population have been confirmed to have contracted COVID-19.",
      "Occurrence is inversely proportional to Rate.",
      "No adjustment has been made in this table for patients who have recovered or passed away.",
      "Whilst the virus spreads with no regard to population size, the rate as shown here indicates the degree of stress which each country's health-care system is likely to be experiencing.",
      "Several countries shown here have not attracted much media attention due lower case numbers than China and Italy, but are likely to be experiencing significant duress.",
      "Continent and Region information is as provided by the UN classifications."
    )
  ) %>%
  formatCurrency(
    columns = c("Population (millions)"),
    currency = "", 
    digits = 3,
    mark = ","
  ) %>%
  formatCurrency(
    columns = c("Confirmed", "Rate (Cases per Million)", "Occurrence"),
    currency = "", 
    digits = 0,
    mark = ","
  )

Cumulative Incidence Rates

ausDays <- confirmed %>% 
  dplyr::filter(Country == "Australia") %>%
  group_by(Country, date) %>%
  summarise_at(vars(confirmed), sum) %>%
  left_join(wikiPops) %>%
  mutate(rate = 1e6 * confirmed / Population) %>%
  dplyr::filter(rate > startingPoint) %>%
  nrow()
minDays <- ausDays - 30
# Use Singapore as that has the longest dataset besides Hubei
nDays <- confirmed %>%
  dplyr::filter(Country == "Singapore") %>% 
  group_by(Country, date) %>%
  summarise(confirmed = sum(confirmed)) %>%
  ungroup() %>%
  left_join(wikiPops) %>%
  mutate(
    rate = 1e6*confirmed / Population
  ) %>% 
  dplyr::filter(rate > startingPoint) %>% 
  nrow() %>%
  subtract(1)

refRate <- c(2, 4, 8)
minPop <- 4e6
p <- confirmed %>%
  group_by(Country, date) %>%
  summarise(confirmed = sum(confirmed)) %>%
  ungroup() %>%
  inner_join(
    dplyr::filter(
      wikiPops, Population > minPop | Country %in% alwaysShow
    )
  ) %>%
  mutate(
    rate = 1e6*confirmed / Population
  ) %>%
  dplyr::filter(
    rate > startingPoint
  ) %>%
  group_by(Country) %>%
  mutate(
    days = date - min(date)
  ) %>%
  dplyr::filter(
    max(days) >= minDays | Country %in% alwaysShow
  ) %>%
  ungroup() %>%
  mutate(
    days = as.integer(days),
    rate = round(rate, 2)
  ) %>%
  dplyr::filter(days <= nDays | Country %in% alwaysShow) %>%
  arrange(date) %>%
  mutate(
    Country = fct_inorder(Country),
    Region = fct_lump(Region, n = 11)
    ) %>%
  rename_all(str_to_title) %>%
  mutate(ymax = max(Rate)) %>%
  ggplot(
    aes(Days, Rate, colour = Country, Date = Date, Confirmed = Confirmed)
  ) +
  geom_segment(
    aes(x, y, xend = xmax, yend = ymax),
    data = . %>%
      dplyr::slice(seq_along(refRate)) %>%
      dplyr::select(ymax) %>%
      mutate(
        ymax = ymax*1.2,
        x = 0,
        y = startingPoint,
        xmax =  refRate*log2(ymax / startingPoint)
      ),
    inherit.aes = FALSE,
    colour = "grey80",
    linetype = 2
  ) +
  geom_line() +
  scale_x_continuous(
    expand = expansion(mult = c(0, 0.05)),
  ) +
  scale_y_log10(
    expand = expansion(mult = c(0, 0.05)),
    label = comma
  ) +
  xlab(
    paste(
      "Days since passing", 
      startingPoint,
      "confirmed cases/million"
    )
  ) +
  ylab("Confirmed Cumulative Infection Rate (cases/million)") +
  facet_wrap(~Region, ncol = 3)

ggplotly(
  p, 
  tooltip = c(
    "Days", "Rate", "Country", "Date", "Confirmed"
  )) 

COVID-19 Confirmed Cumulative Infection Rate for countries which have exceeded 4 confirmed cases/million for 59 or more days, and with populations greater than 4,000,000, apart from a small number of specifically included countries. Data is only shown for the first 122 calendar days since passing 4 confirmed cases/million. Note that from the day records begin in this dataset (2020-01-22), the confirmed infection rate in Hubei was 7.5 confirmed cases/million. Diagonal grey lines indicate a doubling in the infection rate every 2, 4, or 8 days. To hide a country, click on the country in the plot legend. Clicking again on the country in the legend will restore the data within the plot. Countries are shown in order of the date at which they passed the 4 confirmed case/million mark. Regions are as defined by the UN with some regions lumped together if only a few data points were available. Due to the large number of countries shown, you may need to scroll through the legend. Regions of the plot are also able to be zoomed interactively. Please note the y-axis is shown on the logarithmic scale, so that a series of points which appear to be diagonal will indicate exponential growth. The flatter the line, the slower the growth and a perfectly horizontal line would indicate zero growth, or no new confirmed cases.

All figures and tables presented here simply aim to show an alternative viewpoint on the data. Every way to view COVID-19 data will mask important features, and the values shown here do not take into account vital factors such as population density, variability of infection across regions within countries, social culture and demographics. Many countries may not be directly comparable for a combination of the above factors. It is simply to view the data through the lens of a country’s population size using a value which should be easily interpretable.

In the above plot:

• Only countries are shown where the infection rate has surpassed 4 cases/million for more than 59 days. This choice seemed to be a useful way of enabling the comparison of COVID-19 progression across countries after scaling for population size

Daily Increase Vs Confirmed Cases

An alternate viewpoint on the data is to remove time and inspect the relationship between the daily increase in cases and the total number of cases. When this relationship ceases it’s near linear relationship, this can be a sign the control measures have begun to take effect. Whilst this relationship appears to have broken down for Australia, no breakdown has yet occurred for countries such as the Singapore, Canada, the UK, USA and Turkey, with these countries explicitly still in the exponential growth phase.

minPop <- 4e6
minDays <- 15
minRate <- 9
plotIncVConf <- confirmed %>%
  group_by(Country, date) %>%
  summarise_at(vars(confirmed), sum) %>%
  dplyr::filter(confirmed > 0) %>%
  inner_join(
    dplyr::filter(wikiPops, Population > minPop | Country %in% alwaysShow)
  ) %>%
  mutate(
    rate = 1e6 * confirmed / Population
  ) %>%
  split(f = .$Country) %>% 
  lapply(function(x, n = 7){
    x %>%
      mutate(
        d = c(0, diff(rate))
      ) %>%
      mutate_at(
        vars(d), ema, n = n
      ) %>%
      dplyr::filter(
        !is.na(d),
        rate > minRate
      ) %>%
      mutate(
        nDays = nrow(.)
      )
  }
  ) %>%
  bind_rows() %>%
  ungroup() %>%
  dplyr::filter(
    nDays > minDays,
    d > 0.01
  ) %>%
  arrange(date) %>%
  mutate(
    Country = fct_inorder(Country),
    rate = round(rate, 2),
    d = round(d, 3),
    Population = comma(round(Population, -3)),
    Region = fct_lump(Region, n = 11)
  ) %>%
  rename(
    Rate = rate,
    `Average Daily Increase` = d,
    Date = date
  ) %>%
  ggplot(
    aes(Rate, `Average Daily Increase`, colour = Country, label = Population, key = Date)
  ) +
  geom_line() +
  scale_y_log10(
    name = "Average Daily Increase (Cases / Million)",
    label = label_comma(accuracy = 0.1)) +
  scale_x_log10(
    name = "Cumulative Confirmed Infection Rate (Cases / Million)"
  ) +
  facet_wrap(~Region, ncol = 3) 

capIncVConf <- glue(
  "*Daily Increase in Cases plotted against Confirmed Cases, using confirmed cases / million.
  These two values are directly proportional until interventions are successful, at which point the proportional relationship changes, as evidenced by a sudden downwards turn.
  Scaling by population aids in the visualisation of where in the relative infection trajectory each country's control measures have begun to take effect.
  Daily increases are shown using a 7-day exponential moving average in order to minimise the impact of day-to-day variation.
  Countries are only shown from the point at which the moving average exceeds {minRate} cases/million, and have exceeded this value for > {minDays} days.
  Regions are as defined by the UN with some combined for convenience if only a small number of countries were available.
  Data is additionally restricted to countries with a population > {comma(minPop)}.*
  "
)

ggplotly(plotIncVConf +
  coord_cartesian(
    ylim = c(0.1, max(plotIncVConf$data$`Average Daily Increase`)))
)

Daily Increase in Cases plotted against Confirmed Cases, using confirmed cases / million. These two values are directly proportional until interventions are successful, at which point the proportional relationship changes, as evidenced by a sudden downwards turn. Scaling by population aids in the visualisation of where in the relative infection trajectory each country’s control measures have begun to take effect. Daily increases are shown using a 7-day exponential moving average in order to minimise the impact of day-to-day variation. Countries are only shown from the point at which the moving average exceeds 9 cases/million, and have exceeded this value for > 15 days. Regions are as defined by the UN with some combined for convenience if only a small number of countries were available. Data is additionally restricted to countries with a population > 4,000,000.

Dimensional Reduction

In order to summarise which countries are the most similar to each other, a Principal Component Analysis was performed. This enables the multi-dimensional data of the above plots to summarised in two dimensions. As missing data cannot be included in this analysis, several countries which are at earlier comparative time-points than Australia were omitted from this analysis.

nDays <- confirmed %>% 
  dplyr::filter(Country == "Australia") %>% 
  group_by(Country, date) %>%
  summarise(confirmed = sum(confirmed)) %>%
  ungroup() %>%
  inner_join(wikiPops) %>%
  mutate(
    rate = 1e6*confirmed / Population
  ) %>% 
  dplyr::filter(rate > startingPoint) %>% 
  nrow() %>%
  subtract(4)

minPop <- 4e6
pca <- confirmed %>%
  group_by(Country, date) %>%
  summarise(confirmed = sum(confirmed)) %>%
  ungroup() %>%
  inner_join(wikiPops) %>%
  mutate(
    rate = 1e6*confirmed / Population
  ) %>%
  dplyr::filter(
    rate > startingPoint,
    Population > minPop | Country %in% alwaysShow
  ) %>%
  group_by(Country) %>%
  mutate(
    days = date - min(date)
  ) %>%
  dplyr::filter(
    max(days) >= nDays
  ) %>%
  ungroup() %>%
  mutate(
    days = as.integer(days),
    rate = round(rate, 2)
  ) %>%
  dplyr::filter(days <= nDays) %>%
  dplyr::select(Country, rate, days) %>%
  pivot_wider(
    values_from = rate,
    names_from = days
  ) %>%
  as.data.frame() %>%
  column_to_rownames("Country") %>%
  as.matrix() %>%
  .[!rowAnyNAs(.),] %>%
  log() %>%
  prcomp() 

set.seed(101)
pca$x %>%
  as.data.frame() %>%
  rownames_to_column("Country") %>%
  mutate(
    Cluster = cbind(PC1, PC2) %>%
      apply(2, function(x){
        (x - mean(x)) / sd(x)
      }) %>%
      kmeans(centers = k, nstart = 10) %>%
      .[["cluster"]] %>%
      as.factor()
  ) %>%
  ggplot(aes(PC1, PC2, label = Country, colour = Cluster)) +
  geom_point() +
  geom_text_repel(
    show.legend = FALSE
  ) +
  stat_ellipse(
    aes(fill = Cluster),
    colour = NA,
    geom = "polygon",
    alpha = 0.1
    ) +
  xlab(
    paste0(
      "PC1 (", pcPercent[['PC1']],")"
      )
  ) +
    ylab(
    paste0(
      "PC2 (", pcPercent[['PC2']], ")"
      )
  ) +
  theme(
    legend.position = "none"
  )

Dimensional reduction showing which countries infection rates are the most similar to each other at the 85 day time point, after passing 4 confirmed cases/million. This may or may not be indicative of future spread within the population. The value 85 days was simply chosen as this marks how long since the US passed this threshold. Countries which have not progressed beyond 4 confirmed cases/million for 85 days or more are not shown. Countries with populations < 4,000,000 are also excluded. Clusters were generated using k-means, manually specifying k = 3 and are not definitive. Principal Component 1, on the x-axis, corresponds to the greatest source of variability within the data (90.2%), and countries which appear near each other along this axis can be assumed to be showing similar growth in confirmed infection rates at this time point. Separation along the y-axis is less significant, but also worthy of note, as this represents 7.8% of variability within the data.

Fatalities

fr <- confirmed %>%
  inner_join(deaths) %>%
  group_by(Country) %>%
  dplyr::filter(date == max(date)) %>%
  ungroup() %>%
  summarise(fr = sum(deaths) / sum(confirmed)) %>%
  .[["fr"]]

The current fatality rate from all confirmed cases is 5.9%. This may be a function of under-reporting of true cases and is very likely to be an overestimate.

Table of Fatalities

minPop <- 4e6
minCases <- 1000
deaths %>%
  left_join(confirmed) %>%
  group_by(Country, date) %>%
  summarise_at(vars(deaths, confirmed), sum) %>%
  dplyr::filter(deaths > 0) %>%
  left_join(wikiPops) %>%
  ungroup() %>%
  mutate(
    infectionRate = round(1e6 * confirmed / Population, 1),
    fatalityRate = deaths / confirmed,
    Population = round(Population / 1e6, 3)
  ) %>%
  arrange(desc(date), fatalityRate) %>%
  distinct(Country, .keep_all = TRUE) %>%
  arrange(desc(fatalityRate)) %>%
  dplyr::select(
    Continent, Region, Country, Population,
    `Confirmed Cases` = confirmed, 
    `Infection Rate (per million)` = infectionRate,
    `Fatalities` = deaths, 
    `Fatality Rate` = fatalityRate
  ) %>%
  datatable(
    options = list(
      pageLength = 25, 
      autoWidth = TRUE,
      searchCols = list(
        NULL, NULL, NULL, 
        list(
          search = glue(
            '{minPop/1e6} ... {max(wikiPops$Population/1e6)}'
          )
        ),
        list(
          search = glue(
            '{minCases} ... {max(.$`Confirmed Cases`)}'
          )
        ),
        NULL, NULL, NULL
      )
    ),
    filter = 'top',
    class = "stripe",
    rownames = FALSE,
    caption = glue(
      "All countries with recorded fatalities, sorted by default in decreasing order of the Fatality Rate.
      The Fatality Rate simply indicates the number of confirmed cases which end in a fatality.
      The Infection Rate represents the cumulative number of cases confirmed within the country, per million people.
      All columns are searchable and filters are set by default to exclude countries with populations below {comma(minPop)} and those with fewer than {minCases} confirmed cases."
    )
  ) %>%
  formatRound(
    columns = "Infection Rate (per million)",
    mark = ",",
    digits = 1
  ) %>%
  formatPercentage(
    columns = "Fatality Rate",
    digits = 2
  )

Cumulative Fatalities (Unscaled)

minDeaths <- 10
ausDays <- deaths %>%
  dplyr::filter(Country == "Australia") %>% 
  group_by(date) %>%
  summarise_at(vars(deaths), sum) %>% 
  dplyr::filter(deaths > 10) %>%
  nrow()
minDays <- 2
minPop <- 2e6
deathPlot <- confirmed %>%
  left_join(deaths) %>%
  group_by(Country, date) %>%
  summarise_at(
    vars(confirmed, deaths),
    sum
  ) %>%
  ungroup() %>%
  dplyr::filter(
    deaths >= minDeaths,
    Country %in% c(dplyr::filter(wikiPops, Population > minPop)$Country, alwaysShow)
  ) %>%
  group_by(Country) %>%
  mutate(Days = as.integer(date - min(date))) %>%
  dplyr::filter(
    max(Days) > minDays | Country %in% alwaysShow
  ) %>%
  ungroup() %>%
  rename_all(str_to_title) %>%
  arrange(desc(Days)) %>%
  left_join(
    wikiPops %>% dplyr::select(Country, Continent, Region)
  ) %>%
  mutate(
    Country = fct_inorder(Country),
    Region = fct_lump(Region, n = 11)
  ) %>%
  ggplot(
    aes(
      x = Days, y = Deaths, colour = Country,
      label = Date, conf = Confirmed)
  ) +
  geom_line() +
  scale_y_log10(labels = comma) +
  labs(
    x = glue("Days Since Passing {minDeaths} Deaths"),
    y = "Cumulative Fatalities"
  ) + 
  facet_wrap(~Region, ncol = 3)

ggplotly(deathPlot)

As with the spread of the virus, fatalities also grow at an exponential rate. Any slowing in the growth of fatalities is an accurate marker for when the spread of the virus has definitively slowed, despite being a significantly lagging marker. Cumulative fatalities are only shown for countries which have seen 10 or more fatalities for a period beyond 2 days, and for countries with a population > 2,000,000. Regions are shown as defined by the UN, however some with fewer points were merged for convenience.

Cumulative Fatalities (Scaled)

minRate <- 1
minDays <- 24
minPop <- 4e6
scaledDeathPlot <- deaths %>% 
  group_by(Country, date) %>% 
  summarise_at(vars(deaths), sum) %>%
  ungroup() %>%
  dplyr::filter(deaths > 0) %>%
  right_join(
    wikiPops %>%
      dplyr::filter(Population > minPop | Country %in% alwaysShow)
  ) %>%
  mutate(
    rate = round(1e6 * deaths / Population, 2)
  ) %>%
  dplyr::filter(rate > minRate) %>%
  group_by(Country) %>%
  mutate(Days = as.integer(date - min(date))) %>%
  dplyr::filter(
    max(Days) > minDays | Country %in% alwaysShow
  ) %>%
  arrange(desc(Days)) %>%
  ungroup() %>%
  mutate(
    Country = fct_inorder(Country),
    Region = fct_lump(Region, n = 11),
    `One Death Every` = round(Population / deaths, 0),
    deaths = comma(deaths),
    `Population (millions)` = round(Population / 1e6, 1)
  )  %>%
  rename_all(str_to_title) %>%
  rename_all(str_replace_all, pattern = "Mill", replacement = "mill") %>%
  ggplot(
    aes(
      x = Days, y = Rate,  colour = Country, key = Date, 
      b = `Population (millions)`, label = Deaths, 
      family = `One Death Every`
    )
  ) +
  geom_line() +
  scale_y_log10() +
  labs(
    x = glue("Days Since Passing {minRate} Deaths / Million"),
    y = "Deaths / Million"
  ) +
  facet_wrap(~Region, ncol = 3)

scaledDeathPlot %>%
  ggplotly(
    tooltip = c("Country", "Date", "Rate", "Deaths", "One Death Every", "Population (millions)")
  )

Cumulative fatalities for all countries with a population greater than 4 million, who passed 1 deaths / million more than 24 days ago. In the legend at right, countries are shown in order of passing 1 deaths / million. Regions are shown as defined by the UN, however some with fewer points may have been merged for convenience.

Daily Increase in Fatalities Vs Confirmed Fatalities

minCases <- 10
minPop <- 4e6
minDays <- 14
dailyVsTotalDeaths <- deaths %>%
  group_by(Country, date) %>%
  summarise_at(vars(deaths), sum) %>%
  dplyr::filter(deaths > 0) %>%
  split(f = .$Country) %>%
  lapply(function(x, n = 7){
    x %>%
      mutate(
        d = c(0, diff(deaths))
      ) %>%
      mutate_at(
        vars(deaths, d), ema, n = n
      )
  }
  ) %>%
  bind_rows() %>%
  inner_join(wikiPops) %>%
  dplyr::filter(deaths > minCases) %>%
  mutate_at(vars(deaths, d), round, 2) %>%
  mutate(n = n()) %>%
  dplyr::filter(
    n > minDays | Country %in% alwaysShow,
    Population > minPop,
    d > 0
  ) %>%
  rename_all(str_to_title) %>%
  rename(
    `Average Daily Fatalities` = D,
    `Average Total Fatalities` = Deaths
  ) %>%
  arrange(Date) %>%
  ungroup() %>%
  mutate(
    Population = comma(round(Population, -3)),
    Country = fct_inorder(Country),
    Region = fct_lump(Region, n = 11)
    ) %>%
  ggplot(aes(`Average Total Fatalities`, `Average Daily Fatalities`, colour = Country, label = Date, key = Population)) +
  geom_line() +
  scale_x_log10(labels = comma) +
  scale_y_log10(labels = label_comma(accuracy = 0.1)) +
  labs(
    x = "Total Confirmed Fatalities",
    y = 'Daily Fatalities'
  ) +
  facet_wrap(~Region, ncol = 3)

capDailyVsTotal <- glue(
  "*Daily Fatalities plotted against Total Fatalities.
  These two values are directly proportional until interventions are successful, at which point the proportional relationship changes, as evidenced by a sudden downwards turn.
  Values shown are 7-day exponential moving averages in order to minimise the impact of day-to-day variation.
  Countries are only shown from the point at which the moving average exceeds {minCases} cases, and have exceeded this value for > {minDays} days.
  Data is additionally restricted to countries with a population > {comma(minPop)}, for greater clarity.
  Importantly, __due to the time taken from the initial infection to the day of death, this is a lag indicator of the control of infection__.*
  "
)

ggplotly(
  dailyVsTotalDeaths +
    coord_cartesian(
      ylim = c(0.1, max(dailyVsTotalDeaths$data$`Average Daily Fatalities`))
    )
)

Daily Fatalities plotted against Total Fatalities. These two values are directly proportional until interventions are successful, at which point the proportional relationship changes, as evidenced by a sudden downwards turn. Values shown are 7-day exponential moving averages in order to minimise the impact of day-to-day variation. Countries are only shown from the point at which the moving average exceeds 10 cases, and have exceeded this value for > 14 days. Data is additionally restricted to countries with a population > 4,000,000, for greater clarity. Importantly, due to the time taken from the initial infection to the day of death, this is a lag indicator of the control of infection.

Comparison of Total Fatalities

minPop <- 5e6
nCnt <- 25
topDeaths <- deaths %>%
  dplyr::filter(deaths > 0) %>%
  group_by(Country, date) %>%
  summarise_at(vars(deaths), sum) %>%
  dplyr::filter(
    deaths == max(deaths),
    date == max(date)
  ) %>%
  ungroup() %>%
  inner_join(wikiPops) %>%
  # dplyr::filter(
  #   Population > minPop | Country %in% alwaysShow
  # ) %>%
  mutate(
    rate = 1e6*deaths / Population,
    totalRank = nrow(.) - rank(deaths) + 1,
    rateRank = nrow(.) - rank(rate) + 1,
    aveRank = (totalRank + rateRank)/2
    ) %>%
  arrange(aveRank, desc(Population)) %>%
  dplyr::slice(1:nCnt) %>%
  arrange(desc(rateRank)) %>%
  mutate(Country = fct_inorder(Country)) %>%
  arrange(rateRank) %>%
  mutate(Continent = fct_inorder(as.character(Continent)))
a <- topDeaths %>%
  ggplot(aes(Country, rate, fill = Country)) +
  geom_col(colour = "black") +
  geom_label(
    aes(y = rate + 50, label = round(rate, 0)),
    fill = "white", alpha = 0.5
  ) +
  scale_fill_viridis_d(option = "magma", direction = 1) +
  scale_y_continuous(expand = expansion(c(0, 0.08))) +
  labs(y = "Fatalities / Million") +
  coord_flip() +
  facet_grid(Continent~., scales = "free_y", space = "free_y") +
  theme(legend.position = "none")
b <- topDeaths %>%
  ggplot(aes(Country, deaths, fill = Country)) +
  geom_col(colour = "black") +
  geom_label(
    aes(y = deaths + 10000, label = comma(deaths, accuracy = 1)),
    fill = "white", alpha = 0.5
  ) +
  scale_fill_viridis_d(option = "magma", direction = 1) +
  scale_y_continuous(expand = expansion(c(0, 0.08)), labels = comma) +
  labs(y = "Total Fatalities") +
  coord_flip() +
  facet_grid(Continent~., scales= "free_y", space = "free_y") +
  theme(legend.position = "none")
cp <- glue(
  "*Comparison of the {nrow(a$data)} most impacted countries.
  Countries were ranked by total fatalities and by fatalities / million, with the {nrow(a$data)} most highly ranked across both methods are shown.
  Fatalities are shown using A) the number of fatalities scaled by population size (Fatalities / Million) and B) Total Fatalities.
  Countries are also grouped by their continent as designated by the UN classifications, and within each continent each country is ordered by Fatalities / Million.
  For fatality rates scaled by population size, it is important to realise that a value of 500 indicates that one in every 2000 people from the total population has died (ignoring demographics).
  Similarly, a value of 100 indicates that 1 in every 10,000 from the total population has died in that country.
  Whilst the US currently has the most fatalities, 1 in every {comma(round(1e6 / dplyr::filter(a$data, Country == 'US')$rate, 0))} from the US population have currently been confirmed to have died from COVID19.
  For {as.character(a$data$Country[which.max(a$data$rate)])}, 1 in every {comma(round(1e6 / max(a$data$rate), 0))} people are recorded as having died from COVID 19.
  Importantly, whilst the fatality count in countries like the US and the UK are below the statistical value of 'excess deaths' observed throughout these countries, this is not true for Belgium and the fatality rate for Belgium is more likely to reflect an accurate assessment of the true fatalities due to COVID-19.
  In contrast, the values for the US and the UK are likely to be under-estimates of the true fatalities.*"
)

plot_grid(
  a + theme(legend.position = "none"),
  b + 
    theme(
      legend.position = "none",
      axis.title.y = element_blank()
    ),
  labels = c("A", "B"),
  nrow = 1
)

Comparison of the 25 most impacted countries. Countries were ranked by total fatalities and by fatalities / million, with the 25 most highly ranked across both methods are shown. Fatalities are shown using A) the number of fatalities scaled by population size (Fatalities / Million) and B) Total Fatalities. Countries are also grouped by their continent as designated by the UN classifications, and within each continent each country is ordered by Fatalities / Million. For fatality rates scaled by population size, it is important to realise that a value of 500 indicates that one in every 2000 people from the total population has died (ignoring demographics). Similarly, a value of 100 indicates that 1 in every 10,000 from the total population has died in that country. Whilst the US currently has the most fatalities, 1 in every 3,015 from the US population have currently been confirmed to have died from COVID19. For Belgium, 1 in every 1,206 people are recorded as having died from COVID 19. Importantly, whilst the fatality count in countries like the US and the UK are below the statistical value of ‘excess deaths’ observed throughout these countries, this is not true for Belgium and the fatality rate for Belgium is more likely to reflect an accurate assessment of the true fatalities due to COVID-19. In contrast, the values for the US and the UK are likely to be under-estimates of the true fatalities.

Fatality Rate Vs Time

minDays <- 30
minPop <- 10e6
df <- confirmed %>%
  inner_join(deaths) %>%
  group_by(Country, date) %>%
  summarise_at(
    vars(confirmed, deaths),
    sum
  ) %>%
  inner_join(
    wikiPops %>% dplyr::filter(Population > minPop | Country %in% alwaysShow)
  ) %>%
  ungroup() %>%
  mutate(
    `Infection Rate` = 1e6 * confirmed / Population
  ) %>%
  dplyr::filter(
    `Infection Rate` > startingPoint
    ) %>%
  group_by(Country) %>%
  mutate(
    Days = date - min(date),
    n = n()
  ) %>%
  dplyr::filter(
    n > minDays,
    max(deaths, na.rm = TRUE) > 0
  ) %>%
  mutate(
    Rate = deaths / confirmed,
    `Fatality Rate` = percent(Rate, accuracy = 0.1),
    minusT = date - max(date)
  ) %>%
  ungroup() 
plotFr <- mutate(
  df,
  Country = factor(
    Country, 
    levels = df %>%
      dplyr::select(Country, minusT, Rate) %>%
      pivot_wider(
        id_cols = Country, 
        names_from = minusT, 
        values_from = Rate
      ) %>% 
      as.data.frame() %>%
      column_to_rownames("Country") %>% 
      dist() %>% 
      hclust() %>% 
      as.dendrogram() %>% 
      labels()
  )
) %>% 
  rename_all(str_to_title) %>%
  ggplot(
    aes(
      x = Days, y = Country, fill = Rate,
      conf = Confirmed, 
      deaths = Deaths,
      date = Date,
      label = `Fatality Rate`
      )
  ) +
  geom_raster() +
  geom_vline(
    aes(xintercept = Days + 0.5),
    data = . %>%
      dplyr::filter(Country == "Australia") %>%
      dplyr::filter(Date == max(Date)),
    linetype = 3,
    size = 1/3,
    colour = "grey70"
  ) +
  scale_fill_viridis_c(
    option = "magma"
  ) +
  scale_x_continuous(
    expand = expansion(0, 0),
    labels = abs
    ) +
  scale_y_discrete(expand = expansion(0, 0)) +
  labs(
    x = glue(
      "Days since passing {startingPoint} cases/million"
    ),
    y = c(),
    fill = "Fatality\nRate"
  ) +
  theme(
    panel.grid = element_blank()
  )

cpFr <- glue(
  "*Fatality Rate for confirmed cases after passing {startingPoint} confirmed cases / million.
  Only countries with {minDays} days of data beyond this time-point and a population size >{comma(minPop)} are shown.
  Country order is based on clustering using the most recent values.
  Countries with no recorded fatalities have been excluded for obvious reasons. 
  The dashed grey line indicates the time-point Australia is currently at.
  A clear trend of an increasing fatality rate with time is evident in many countries (e.g. Spain, France, Italy, UK), however, for some countries (e.g. Singapore) this rate appears relatively stable throughout the majority of days.
  The overall Fatality Rate for confirmed cases is currently ({percent(fr, accuracy = 0.1)}).*"
)

ggplotly(
  plotFr ,
  tooltip = c(
    "Country", "Date", "Days", "Confirmed", "Deaths",
    "Fatality Rate"
    )
  )

Fatality Rate for confirmed cases after passing 4 confirmed cases / million. Only countries with 30 days of data beyond this time-point and a population size >10,000,000 are shown. Country order is based on clustering using the most recent values. Countries with no recorded fatalities have been excluded for obvious reasons. The dashed grey line indicates the time-point Australia is currently at. A clear trend of an increasing fatality rate with time is evident in many countries (e.g. Spain, France, Italy, UK), however, for some countries (e.g. Singapore) this rate appears relatively stable throughout the majority of days. The overall Fatality Rate for confirmed cases is currently (5.9%).

Currently Active Infections

rr <- confirmed %>%
  group_by(Country, date) %>%
  summarise_at(vars(confirmed), sum, na.rm = TRUE) %>%
  left_join(
    recovered %>% 
      group_by(Country, date) %>%
      summarise_at(vars(recovered), sum, na.rm = TRUE)
  ) %>%
  dplyr::filter(date == max(date)) %>%
  ungroup() %>%
  summarise(rr = sum(recovered) / sum(confirmed)) %>%
  .[["rr"]]

Information regarding recovered cases is likely to be the least reliable of reported values as many regions do not report updated numbers for several consecutive days. Additionally many regions do not report recovered cases as the criteria for considering a person to have recovered as currently unclear. Given this:

• The current recovery rate from all confirmed cases is 40.8%. Given the level of under-reporting this may also be highly inaccurate.
• At this point, 53.3% of all confirmed infections are considered as ‘active’.

Active Infections Across All Countries

minDays <- 25
minPop <- 4e6
df <- confirmed %>%
  dplyr::filter(
    confirmed > 0,
    Country != "China (Other)"
  ) %>%
  left_join(deaths) %>%
  group_by(Country, date) %>%
  summarise_at(
    vars(confirmed, deaths), sum
  ) %>%
  left_join(
    recovered %>%
      group_by(Country, date) %>%
      summarise_at(vars(recovered), sum)
  ) %>%
  mutate_at(vars(confirmed, deaths, recovered), cummax) %>%
  mutate(
    active = confirmed - deaths - recovered
  ) %>%
  dplyr::filter(!is.na(active)) %>%
  inner_join(
    dplyr::filter(wikiPops, Population > minPop | Country %in% alwaysShow)
  ) %>%
  mutate(rate = 1e6 * active / Population) %>%
  dplyr::filter(rate > startingPoint) %>%
  group_by(Country) %>%
  mutate(
    days = date - min(date)
  ) %>%
  dplyr::filter(max(days) > minDays | Country %in% alwaysShow) %>% 
  ungroup() %>%
  mutate(
    days = as.integer(days),
    rate = round(rate, 2)
  ) %>%
  arrange(date) %>%
  mutate(
    Country = fct_inorder(Country),
    Region = fct_lump(Region, n = 11)
  ) %>%
  rename_all(str_to_title)
nDays <- max(df$Days)
p2 <- df %>%
  ggplot(
    aes(
      x = Days, y = Rate, colour = Country, 
      Date = Date, Active = Active,
      Confirmed = Confirmed, Recovered = Recovered,
      Deaths = Deaths
      )
  ) +
  geom_segment(
    aes(x, y, xend = xmax, yend = ymax),
    data = tibble(
      x = 0,
      y = startingPoint,
      xmax = ceiling(nDays / c(4, 2, 1)),
      ymax = startingPoint + 2^(xmax/ c(2, 4, 8))
    ),
    inherit.aes = FALSE,
    colour = "grey90",
    linetype = 2
  ) + 
  geom_line() +
  scale_x_continuous(
    expand = expansion(mult = c(0, 0.05)),
  ) +
  scale_y_log10(
    expand = expansion(mult = c(0, 0.05))
  ) +
  xlab(
    paste(
      "Days since passing", 
      startingPoint, 
      "confirmed active cases/million"
    )
  ) +
  ylab("Confirmed Active Infection Rate (cases/million)") + 
  facet_wrap(~Region, ncol = 3)

ggplotly(
  p2 +
    coord_cartesian(ylim = range(p2$data$Rate))
) 

Confirmed active cases of COVID-19 for countries where the confirmed infection rate has exceeded 4 confirmed active cases/million for more than 122 calendar days. Only countries with a population greater than 4,000,000 are shown for better visualisation. Due to difficulties introduced by the currently reported low active infection rate outside Hubei province, data from China has been excluded from this plot, with the exception of Hubei and Hong Kong. Recovered cases as poorly and irregularly reported by many countries and as such, this plot may be suffer multiple instances of a sudden decline. To hide a country, click on the country in the plot legend. Clicking again on the country in the legend will restore the data within the plot. Countries are shown in order of the date at which they passed the 4 confirmed active case/million mark. Due to the number of countries shown, you may need to scroll through the legend. Regions of the plot are also able to be zoomed interactively. Please note the y-axis is shown on the logarithmic scale, so that a series of points which appear to be diagonal will indicate exponential growth/decay. Given the different starting point to the previous plot, data will generally be shown for fewer time-points.

Notable features of this perspective are:

• Active infection rates in many countries are currently declining especially in Australia and New Zealand
• Despite strict control measures early on in the pandemic, active infection rates are still exponentially increasing in Singapore

Summary of All Rates

minPop <- 10e6
p4 <- confirmed %>%
  dplyr::filter(
    confirmed > 0
  ) %>%
  left_join(deaths) %>%
  group_by(Country, date) %>%
  summarise_at(vars(confirmed, deaths), sum) %>%
  dplyr::filter(date == max(date)) %>%
  left_join(
    recovered %>%
      group_by(Country, date) %>%
      summarise_at(vars(recovered), sum)
  ) %>%
  ungroup() %>%
  inner_join(
    wikiPops %>% dplyr::filter(Population > minPop | Country %in% alwaysShow)
  ) %>%
  mutate(rate = 1e6 * confirmed / Population) %>%
  dplyr::filter(rate > startingPoint) %>% 
  group_by(Country) %>%
  mutate(
    active = confirmed - recovered - deaths
  ) %>%
  mutate(
    active = 100*active / confirmed,
    recovered = 100*recovered / confirmed,
    fatalities = 100*deaths / confirmed
  ) %>%
  ungroup() %>%
  dplyr::filter(active < 100) %>%
  arrange(active) %>% 
  mutate(Country = fct_inorder(Country)) %>%
  pivot_longer(
    cols = c(active, recovered, fatalities),
    names_to = "Status",
    values_to = "Percentage"
  ) %>%
  mutate(
    Status = str_to_title(Status),
    Status = factor(
      Status, 
      levels = c("Active", "Recovered", "Fatalities")
    ),
    Percentage = round(Percentage, 2)
  ) %>%
  mutate(confirmed = comma(confirmed)) %>%
  rename(Confirmed = confirmed) %>%
  ggplot(
    aes(
      Country, Percentage,
      fill = Status, cases = Confirmed
    )
  ) +
  geom_col() +
  scale_fill_manual(
    values = c(
      Active = "blue",
      Recovered = "green",
      Fatalities = "red"
      )
  ) +
  scale_y_continuous(expand = expansion(0, 0)) +
  coord_flip() +
  labs(x = c()) +
  theme(
    legend.position = "none"
  )

ggplotly(p4)

Fatality, Recovery and Active Infection rates for countries which have exceeded 4 confirmed cases / million, and with a population size > 10,000,000. Countries are ordered by the percentageof cases that remain active.

Australian Specific Data

ausPops <- tribble(
  ~State, ~Population,
  "New South Wales",    8117976,
  "Victoria", 6629870,
   "Queensland", 5115451,
  "South Australia", 1756494,
  "Western Australia", 2630557,
  "Tasmania", 535500,
  "Northern Territory", 245562,
  "Australian Capital Territory", 428060
)

Australian State populations were taken from the ABS Website and were accurate in Sept 2019. The difference with previous estimates used above was within 0.04%, and as such no adjustments were made.

A series of complimentary charts regarding the spread of COVID-19 are available from the ABC website.

• It is now clear that Australia is no longer in an exponential growth phase and they pandemic may be under control
• Using an estimated population size of 25,203,198, the total percentage of the Australian population confirmed as having been infected currently sits at 0.029%, or one person in every 3,473.
• Notably, the proportion of recovered cases is now >80%

confirmed %>% 
  dplyr::filter(
    Country == "Australia", 
    date >= sort(unique(date), decreasing = TRUE)[2]
  ) %>%
  bind_rows(
    group_by(., date) %>%
      summarise(
        confirmed = sum(confirmed)
      ) %>%
      mutate(
        `Province/State` = "Total"
      )
  ) %>%
  group_by(`Province/State`) %>%
  mutate(
    Increase = c(NA, diff(confirmed)),
    `% Increase` = c(NA, diff(confirmed)) / min(confirmed)
  ) %>%
  ungroup() %>%
  pivot_wider(
    id_cols = `Province/State`, 
    names_from = date, 
    values_from = c(confirmed, Increase, `% Increase`)
  ) %>%
  dplyr::select_if(function(x){sum(is.na(x)) <= 1}) %>%
  rename(State = `Province/State`) %>%
  rename_all(str_remove_all, pattern = "confirmed_") %>%
  rename_all(str_remove_all, pattern = "_2020.+") %>%
  left_join(
    dplyr::select(latestAU, State, confirmed, deaths, recovered)
  ) %>%
  mutate_at(
    vars(confirmed, deaths, recovered),
    function(x){
      x[is.na(x)] <- sum(x, na.rm = TRUE)
      x
    }
  ) %>%
  mutate(
    `% Increase` = percent(`% Increase`, accuracy = 0.01),
    `Fatality Rate` = percent(deaths / confirmed, accuracy = 0.1),
    `Recovery Rate` = percent(recovered / confirmed, accuracy = 0.1),
    `Currently Active` = confirmed - recovered - deaths,
  ) %>%
  rename(
    Fatalities = deaths,
    Recovered = recovered
  ) %>%
  dplyr::select(
    State, starts_with("20"), ends_with("Increase"), starts_with("Fatal"), starts_with("Recov"), `Currently Active`
  ) %>%
  pander(
    justify = "lrrrrrrrrr",
    caption = paste(
      "*Confirmed cases, fatalities and recoveries reported by each state at time of preparation.",
      "Any states with unchanged, or decreasing confirmed cases may indicate delays with the automated data sources, such as health.gov.au or JHU, or that these states have not yet reported for the day.", 
      "Please note that some discrepancy with dates may occur due to automated data sources obtained different time zones, such as the USA, the UK and Australia.*"
    ),
    emphasize.strong.rows = nrow(.)
  )

Confirmed cases, fatalities and recoveries reported by each state at time of preparation. Any states with unchanged, or decreasing confirmed cases may indicate delays with the automated data sources, such as health.gov.au or JHU, or that these states have not yet reported for the day. Please note that some discrepancy with dates may occur due to automated data sources obtained different time zones, such as the USA, the UK and Australia.

State	2020-06-05	2020-06-06	Increase	% Increase	Fatalities	Fatality Rate	Recovered	Recovery Rate	Currently Active
Australian Capital Territory	107	107	0	0.00%	3	2.8%	104	97.2%	0
New South Wales	3,110	3,110	0	0.00%	50	1.6%	2,719	87.4%	341
Northern Territory	30	30	0	0.00%	0	0.0%	30	100.0%	0
Queensland	1,061	1,061	0	0.00%	6	0.6%	1,049	98.9%	6
South Australia	440	440	0	0.00%	4	0.9%	436	99.1%	0
Tasmania	228	228	0	0.00%	13	5.8%	210	92.9%	3
Victoria	1,681	1,681	0	0.00%	19	1.1%	1,586	94.3%	76
Western Australia	596	599	3	0.50%	9	1.5%	559	93.3%	31
Total	7,253	7,256	3	0.04%	104	1.4%	6,693	92.3%	457

Plot of Current Australian Values

ausStatsPlot <- ggplotly(
  confirmed %>%
    dplyr::filter(Country == "Australia", confirmed > 0) %>%
    left_join(deaths) %>%
    left_join(recovered) %>%
    group_by(Country, date) %>%
    summarise_at(vars(confirmed, deaths, recovered), sum) %>%
    ungroup() %>%
    mutate_at(vars(confirmed, deaths, recovered), cummax) %>%
    mutate(active = confirmed - deaths - recovered) %>%
    pivot_longer(
      cols = c(active, confirmed, deaths, recovered),
      names_to = "Type",
      values_to = "Total"
    ) %>%
    arrange(Type, date) %>%
    dplyr::filter(date > ymd("2020-02-29")) %>%
    mutate(
      Type = str_to_title(Type),
      Type = str_replace_all(Type, "Deaths", "Fatalities"),
    ) %>%
    dplyr::filter(Total > 0) %>%
    rename_all(str_to_title) %>%
    ggplot(aes(Date, Total, colour = Type)) +
    geom_point() +
    geom_smooth(method = "loess", se = FALSE, size = 1/2) +
    scale_y_log10() +
    scale_colour_manual(
      values = c(
        Active = rgb(0, 0, 0),
        Confirmed = rgb(0, 0.3, 0.7),
        Fatalities = rgb(0.8, 0.2, 0.2),
        Recovered = rgb(0.2, 0.7, 0.4)
      )
    )
)
ausStatsPlot$x$data[5:8] %<>%
  lapply(function(x){
    x$hoverinfo <- "none"
    x
  })

ausStatsCap <- "*Current confirmed and recovered cases, along with fatalities for Australia only. Active cases are shown as confirmed cases excluding fatalities and those classed as recovered. Loess curves through all points are shown as continuous lines. Data is only shown from 1^st^ March 2020 as this was the date of the first recorded fatality in Australia. Recovered patient information was also sparse in the early stages of data collection, and as a result estimates of active infections will be a significant underestimate until 6^th^ April. In particular, QLD only began reporting recovered cases on this date. NSW followed a fortnight after this date and as such, only the most recent numbers can be considered as accurate. Below this plot, the same figures can be seen broken down by state.*"

ausStatsPlot 

Current confirmed and recovered cases, along with fatalities for Australia only. Active cases are shown as confirmed cases excluding fatalities and those classed as recovered. Loess curves through all points are shown as continuous lines. Data is only shown from 1^st March 2020 as this was the date of the first recorded fatality in Australia. Recovered patient information was also sparse in the early stages of data collection, and as a result estimates of active infections will be a significant underestimate until 6^th April. In particular, QLD only began reporting recovered cases on this date. NSW followed a fortnight after this date and as such, only the most recent numbers can be considered as accurate. Below this plot, the same figures can be seen broken down by state.

ggplotly(
  confirmed %>%
    dplyr::filter(Country == "Australia") %>%
    left_join(deaths) %>%
    left_join(recovered) %>%
    rename(State = `Province/State`) %>%
    dplyr::filter(
      date > ymd("2020-03-19"),
    ) %>%
    arrange(date) %>%
    group_by(State) %>%
    mutate_at(
      vars(confirmed, recovered, deaths), cummax
    ) %>%
    ungroup() %>%
    left_join(ausPops) %>%
    mutate(active = confirmed - deaths - recovered) %>%
    pivot_longer(
      cols = c(confirmed, deaths, recovered, active),
      names_to = "status",
      values_to = "count"
    ) %>%
    dplyr::filter(
      count > 0,
      !(State %in% c("Queensland", "New South Wales") & status == "recovered" & date < ymd("2020-04-06")),
      !(State %in% c("South Australia") & status == "recovered" & date < ymd("2020-04-01")),    
      !(State %in% c("Tasmania") & status == "recovered" & date < ymd("2020-04-02")),
    ) %>%
    mutate(
      rate = 1e6*count/Population,
      rate = round(rate, 2),
      status = str_replace(status, "deaths", "fatal") %>% str_to_title(),
      status = factor(status, levels = c("Confirmed", "Active", "Recovered", "Fatal"))
    ) %>%
    rename_all(str_to_title) %>%
    ggplot(aes(Date, Rate, colour = Status, label = Count)) +
    geom_point() +
    geom_smooth(method = "loess", se = FALSE, size = 1/2) +
    scale_y_log10() +
    scale_colour_manual(
      values = c(
        Active = rgb(0, 0, 0),
        Confirmed = rgb(0, 0.3, 0.7),
        Fatal = rgb(0.8, 0.2, 0.2),
        Recovered = rgb(0.2, 0.7, 0.4)
      )
    ) +
    facet_wrap(~State, ncol = 4) +
    labs(
      x = "Date",
      y = "Cases / Million",
      colour = "Case Status"
    )
)

Breakdown of individual states. Victorian recovered numbers began to be accurately reported from 22^nd March, with other states gradually providing this information. NSW/QLD recovered cases have only recently begun being reported and up until the most recent dates, recovered/active values were very approximate for these states.

State-wise Comparison of Confirmed Infection Rates

minRate <- 10
p5 <- confirmed %>%
  dplyr::filter(
    Country == "Australia"
  ) %>%
  dplyr::rename(State = `Province/State`) %>%
  left_join(ausPops) %>%
  mutate(
    Rate = round(1e6*confirmed / Population, 2),
    Date = format.Date(date, "%d-%B")
  ) %>%
  dplyr::filter(
    !is.na(Population),
    Rate > minRate
  ) %>%
  arrange(date) %>%
  mutate(
    State = fct_inorder(State)
  ) %>%
  dplyr::rename(
    Confirmed = confirmed
  ) %>%
  ggplot(
    aes(
      x = date, y = Rate, colour = State, 
      label = Date, key = Confirmed
      )
  ) +
  geom_point() +
  # geom_smooth(
  #   data = . %>%
  #     dplyr::filter(date >= dt - 14),
  #   method = "lm",
  #   se = FALSE,
  #   show.legend = FALSE,
  #   size = 1/2
  # ) +
  geom_smooth(
    method = "loess",
    se = FALSE,
    linetype = 3,
    size = 1/2
  ) +
  scale_y_log10() +
  labs(
    x = "Date",
    y = "Confirmed Infection Rate (cases/million)"
  )
# Hide the tooltip from the regression lines
n <- length(levels(p5$data$State))
p5 <- ggplotly(p5, tooltip = c("Date", "Rate", "State", "Confirmed"))
regIndex <- seq(n + 1, length.out = length(p5$x$data) , by = 1)
p5$x$data[regIndex] <- lapply(
  p5$x$data[regIndex],
  function(x){
  x$hoverinfo <- "none"
  x
})

p5

Infection rates for each state with data beginning for each state once 10 confirmed cases /million was exceeded. States are shown in order of the initial date they passed 10 cases/million. Most states show a strong upward curve during the week following 19^th March, clearly showing the impact of the Ruby Princess, however the consistent deviation below regression lines for all states is strongly suggestive of a slowing in the recent infection rate. The hypothetical growth in the infection rate without that event occurring is not calculated here, however it is clear that the spread of COVID-19 has not continued at the same pace after that event initially occurred.

Assessment of Testing Within Each State

latestAU %>%
  dplyr::select(State, Tested = tested) %>%
  write_tsv(
    glue(
      "tested/tested_{format(dt, '%Y%m%d')}.tsv"
    )
  )

Testing numbers were initially sourced from manual inspection of individual press releases for NSW, QLD, VIC, WA, SA, TAS, ACT and the NT. Updates on testing numbers beyond the initial values were performed automatically using the above code which probes each state’s official releases for the latest values, and updates where found. The number of tested individuals in each state was then assessed as a function of State population size. All results are valid at the time of report generation.

latestAU %>%
  dplyr::select(State, Tested = tested) %>%
  full_join(
    confirmed %>%
      dplyr::filter(
        date == max(date),
        Country == "Australia"
      ) %>%
      rename(
        State = `Province/State`
      ) 
  ) %>%
  mutate(
    Tested = case_when(
      is.na(Tested) ~ confirmed,
      Tested < confirmed ~ confirmed,
      !is.na(Tested) ~ Tested
    )
  ) %>%
  left_join(ausPops) %>%
  bind_rows(
    tibble(
      State = "Total",
      Population = sum(.$Population, na.rm = TRUE),
      confirmed = sum(.$confirmed, na.rm = TRUE),
      Tested = sum(.$Tested, na.rm = TRUE)
    )
  ) %>%
  mutate(
    `Tests / '000` = round(1e3 * Tested / Population, 2),
    Positive = confirmed / Tested,
    Negative = 1 - Positive,
    isTotal = grepl("Total", State)
  ) %>%
  dplyr::select(
    State, Population,
    Confirmed = confirmed,
    Tested, 
    contains("000"), 
    ends_with("ive"),
    isTotal
  ) %>%
  arrange(isTotal, desc(`Tests / '000`)) %>%
  dplyr::select(-isTotal) %>%
  dplyr::rename(
    `% Positive Tests` = Positive,
    `% Negative Tests` = Negative
  ) %>%
  mutate_at(
    vars(starts_with("%")), percent, accuracy = 0.01
  ) %>%
  pander(
    justify = "lrrrrrr",
    missing = "",
    caption = glue(
      "*COVID-19 testing scaled by state population size.
      Confirmed cases are assumed to be the tests returning a positive result.
      The current numbers available for some states are a lower limit, and as such, the proportion of the population tested is likely to be higher, as is the proportion of tests returning a negative result.*"
    ),
    emphasize.strong.rows = nrow(.)
  )

COVID-19 testing scaled by state population size. Confirmed cases are assumed to be the tests returning a positive result. The current numbers available for some states are a lower limit, and as such, the proportion of the population tested is likely to be higher, as is the proportion of tests returning a negative result.

State	Population	Confirmed	Tested	Tests / ’000	% Positive Tests	% Negative Tests
Victoria	6,629,870	1,681	542,000	81.75	0.31%	99.69%
New South Wales	8,117,976	3,110	555,956	68.48	0.56%	99.44%
Tasmania	535,500	228	35,669	66.61	0.64%	99.36%
South Australia	1,756,494	440	110,000	62.62	0.40%	99.60%
Australian Capital Territory	428,060	107	20,132	47.03	0.53%	99.47%
Western Australia	2,630,557	599	112,110	42.62	0.53%	99.47%
Queensland	5,115,451	1,061	216,458	42.31	0.49%	99.51%
Northern Territory	245,562	30	8,983	36.58	0.33%	99.67%
Total	25,459,470	7,256	1,601,308	62.9	0.45%	99.55%

Current Growth Factor

list(
  confirmed %>%
    dplyr::filter(Country == "Australia") %>%
    arrange(date) %>%
    group_by(
      `Province/State`
    ) %>%
    mutate(
      new = c(0, diff(confirmed)),
      new_ma = ema(new, 7)
    ) %>%
    dplyr::filter(confirmed > 0, !is.na(new_ma)) %>%
    mutate(
      R = c(NA, new_ma[-1] / new_ma[-n()]),
      R = case_when(
        is.nan(R) ~ 0,
        !is.nan(R) ~ R
      )
    ) %>%
    ungroup() %>%
    arrange(`Province/State`),
  confirmed %>%
    dplyr::filter(Country == "Australia") %>%
    arrange(date) %>%
    group_by(Country, date) %>%
    summarise_at(vars(confirmed), sum) %>%
    ungroup() %>%
    mutate(
      new = c(0, diff(confirmed)),
      new_ma = ema(new, 7)
    ) %>%
    dplyr::filter(confirmed > 0, !is.na(new_ma)) %>%
    mutate(
      R = c(NA, new_ma[-1] / new_ma[-n()]),
      R = case_when(
        is.nan(R) ~ 0,
        !is.nan(R) ~ R
      ),
      `Province/State` = "All States"
    ) %>%
    arrange(`Province/State`)
) %>%
  bind_rows() %>%
  dplyr::filter(date > ymd("2020-03-15")) %>%
  ggplot(aes(date, R, colour = `Province/State`)) +
  geom_ribbon(aes(ymin = 1, ymax = R), alpha = 0.1) +
  geom_hline(yintercept = 1) +
  geom_smooth(se = FALSE, linetype = 2, size = 1/3) +
  geom_label(
    aes(label = R),
    data = . %>%
      dplyr::filter(date == max(date)) %>%
      mutate(R = round(R, 2), date = date + 1),
    fill = rgb(1, 1, 1, 0.3),
    show.legend = FALSE
  ) +
  labs(
    x = "Date", y = "Growth Factor"
  ) +
  facet_wrap(~`Province/State`, scales = "free_x") +
  theme(legend.position = "none") +
  coord_cartesian(ylim= c(0, 5))

Growth factor for each State/Territory. This value becomes volatile when daily new cases approach zero as is commonly observed in small populations, and at the end stages of an outbreak. In order to try and minimise volatility a 7 day (i.e. 1 week) exponential moving average was used, in contrast to the simple 5-day average advocated here. Dashed lines are a loess curve

Scaled Transmissions

offset <- 6

An alternative viewpoint to the growth factor is to simply assess the number of new cases arising from the number of active cases. As there is a delay in symptom onset, the number of new cases on a given day is a function of the number of active cases several days in the past, and using an offset value easily allows this to communicated in an informative way. Any sustained number below 1 indicates a trajectory which leads to control of the outbreak, whilst values > 1 indicate exponential growth. This visualisation fails to account for imported cases and assumes that all cases were transmitted within each state. Recent patterns in both WA and SA have shown many cases are still being imported from overseas.

confirmed %>%
  dplyr::filter(Country == "Australia") %>%
  left_join(recovered) %>%
  left_join(deaths) %>%
  mutate(active = confirmed - recovered - deaths) %>%
  bind_rows(
    group_by(., Country, date) %>%
      summarise_at(
        vars(confirmed, recovered, deaths, active),
        sum
      ) %>%
      ungroup() %>%
      mutate(
        `Province/State` = "All States"
      )
  ) %>%
  group_by(`Province/State`) %>%
  mutate(
    new = c(0, diff(confirmed)),
    active_offset = c(rep(NA, offset - 1), active[seq_len(n() - offset + 1)]),
    trans = new / active_offset
  ) %>%
  ungroup() %>%
  dplyr::filter(
    date > ymd("2020-03-09"),
    !is.na(trans)
  ) %>%
  ggplot(aes(date - offset, trans, colour = `Province/State`)) +
  geom_point(aes(x = ymd("2020-04-01"), y = 2.1), colour = NA) +
  geom_vline(
    xintercept = ymd("2020-03-23"), 
    linetype = 2,
    colour = "red"
  ) +
  geom_vline(
    xintercept = ymd("2020-03-15"), 
    linetype = 2,
    colour = "blue"
  ) +
  geom_ribbon(aes(ymax = 1, ymin = trans), alpha = 0.2) +
  geom_hline(yintercept = 1, colour = "grey20") +
  facet_wrap(~`Province/State`, scales = "free_y") +
  labs(
    x = "Transmission Date",
    y = "New Cases / Active Case"
  ) +
  theme(legend.position = "none")

Numbers of new cases arising from existing active cases. The value for active cases was taken using an offset of 6 days to represent the median incubation period and the approximate transmission date is shown on the x-axis. The introduction of social distancing is indicated in blue at the date these effects will have become visible in the data (i.e. 6 days after introduction). The date that the effects of lockdown should have begun appearing is shown in red. Importantly, recovered cases were poorly reported in NSW and QLD until recently, and as such, active cases from early March to early April are a significant overestimate and the ratio of new cases scaled by active cases makes the initial outbreak appear far more moderate than the reality was. Even taking this into account for the actual values, it is clear from the peaks above that social distancing, whilst having a positive impact, was insufficient to fully contain the outbreak.

R Session Information

R version 3.6.3 (2020-02-29)

Platform: x86_64-pc-linux-gnu (64-bit)

locale: LC_CTYPE=C, LC_NUMERIC=C, LC_TIME=C, LC_COLLATE=C, LC_MONETARY=C, LC_MESSAGES=en_AU.UTF-8, LC_PAPER=en_AU.UTF-8, LC_NAME=C, LC_ADDRESS=C, LC_TELEPHONE=C, LC_MEASUREMENT=en_AU.UTF-8 and LC_IDENTIFICATION=C

attached base packages: stats, graphics, grDevices, utils, datasets, methods and base

other attached packages: QuantTools(v.0.5.7), data.table(v.1.12.8), cowplot(v.1.0.0), plotly(v.4.9.2.1), DT(v.0.13), pander(v.0.6.3), RCurl(v.1.98-1.2), rvest(v.0.3.5), xml2(v.1.3.2), jsonlite(v.1.6.1), glue(v.1.4.1), broom(v.0.5.6), ggrepel(v.0.8.2), matrixStats(v.0.56.0), scales(v.1.1.1), lubridate(v.1.7.8), magrittr(v.1.5), forcats(v.0.5.0), stringr(v.1.4.0), dplyr(v.0.8.5), purrr(v.0.3.4), readr(v.1.3.1), tidyr(v.1.0.3), tibble(v.3.0.1), ggplot2(v.3.3.0) and tidyverse(v.1.3.0)

loaded via a namespace (and not attached): httr(v.1.4.1), splines(v.3.6.3), viridisLite(v.0.3.0), here(v.0.1), modelr(v.0.1.7), fasttime(v.1.0-2), assertthat(v.0.2.1), highr(v.0.8), selectr(v.0.4-2), cellranger(v.1.1.0), yaml(v.2.2.1), pillar(v.1.4.4), backports(v.1.1.7), lattice(v.0.20-41), digest(v.0.6.25), colorspace(v.1.4-1), Matrix(v.1.2-18), htmltools(v.0.4.0), pkgconfig(v.2.0.3), haven(v.2.2.0), mgcv(v.1.8-31), generics(v.0.0.2), farver(v.2.0.3), ellipsis(v.0.3.0), withr(v.2.2.0), lazyeval(v.0.2.2), cli(v.2.0.2), crayon(v.1.3.4), readxl(v.1.3.1), evaluate(v.0.14), fs(v.1.4.1), fansi(v.0.4.1), MASS(v.7.3-51.6), nlme(v.3.1-147), Cairo(v.1.5-12), tools(v.3.6.3), hms(v.0.5.3), lifecycle(v.0.2.0), munsell(v.0.5.0), reprex(v.0.3.0), compiler(v.3.6.3), rlang(v.0.4.6), grid(v.3.6.3), rstudioapi(v.0.11), htmlwidgets(v.1.5.1), crosstalk(v.1.1.0.1), bitops(v.1.0-6), labeling(v.0.3), rmarkdown(v.2.1), gtable(v.0.3.0), DBI(v.1.1.0), curl(v.4.3), R6(v.2.4.1), knitr(v.1.28), rprojroot(v.1.3-2), stringi(v.1.4.6), Rcpp(v.1.0.4.6), vctrs(v.0.3.0), dbplyr(v.1.4.3), tidyselect(v.1.1.0) and xfun(v.0.13)